// key ideas:
// offsets are the distances from the beginnings of the current row, page, book (just like a pointer is the distance from beginning of memory)
// step is the distance to get from the end of a row to begin of next row (pitch - flat_size)
// (pitch here is an image processing jargon, it the image size including padding)
// flat_size is the total number of elements in given row/page/book, just like offset
// flat_size can be obtained with (end - begin) of a given range, the relevant calculation(width*height*depth*...) performed when forming the end pointer, unless it might be otherwise formed, through iteration for example
// done_mask (it == end) indicates which dimensions are at their end and need to be stepped forward,
// for example in 2D case
// (false, false) - we are safely inside the range
// (true, false) - we reached the end of the row
// (true, true) - we reached the end of the whole range
// (false, true) - not possible, invariant broken
// in boolean context it reduces by conjunction
// iterator::equality is an optimization of this logic, exploiting the patterns of iteration
// begin != end is similar and thus the following pattern
// while (begin != end) ++begin;
// will only take you till the end of the row,
// afterwords the step will need to be added (next, advance) according to done_mask
// in iterator this is also optimized to a plain bool, so the dimensional information is lost

// TODO: compare to Boost.GIL which stores a pointer and pitch, various advantages and disadvantages, reimplement all the relevant examples from there!!!

#ifndef SIMPLE_GEOM_ITERATOR_HPP
#define SIMPLE_GEOM_ITERATOR_HPP

#include "vector.hpp"

namespace simple::geom
{

	//TODO: in image processing step/pitch can be in bytes not in sizeof(value_type)s, for memory alignment purposes,  which means that offsets will need to be in bytes as well, but maybe that's out of scope here, in general step/pitch is necessary for sub-ranges/regions

	// a more sophisticated and somewhat optimized version (see vectirator)
	// the underlying iterator FlatIt serves as the last offset
	// so the logical structure is
	// 2D - vector(row_offset, flat_iterator)
	// 3D - vector(row_offset, page_offset, flat_iterator)
	// 4D - vector(row_offset, page_offset, book_offset, flat_iterator)
	// etc
	// advantages:
	// performs better in practice,
	// disadvantages:
	// complex implementation
	// very specific interface
	template <typename FlatIt, size_t Dimensions = 2>
	class iterator
	{

		public:

		using offset_type = vector<std::uintptr_t, Dimensions-1>;
		using scalar_difference_type = typename std::iterator_traits<FlatIt>::difference_type;
		using difference_type = vector<scalar_difference_type, Dimensions>;
		using value_type = typename std::iterator_traits<FlatIt>::value_type;
		using reference = typename std::iterator_traits<FlatIt>::reference;
		using pointer = typename std::iterator_traits<FlatIt>::pointer;

		// NOTE: if FlatIt is random access, we are partially random access - adding and subtracting offsets is cheap, but offsets are multidimensional, so the required equivalence with repeated increments doesn't work, that is you can decrement the offset while it's not equal (by conjunction) to default constructed offset, incrementing the pointer, and it will be a valid operation, but not the same as adding the offset directly, the (conjunctive) inequality check will hit a wall on the end of first dimension, while the offset addition can get you across
		// in more general terms something like binary search wouldn't work, cause halfing the offset is meaningless, when pitch is not zero
		// TODO: we are actually input if FlatIt is input, output if output, forward if forward, bidirectional if anything higher, and otherwise not an iterator at all... need the meta-programming here and the proper free functions used on FlatIt member instead of operators
		using iterator_category = std::bidirectional_iterator_tag;
		//TODO: iterator_concept

		constexpr iterator(FlatIt flat, offset_type offsets)
			: offsets(offsets), flat(flat)
		{}
		constexpr iterator(FlatIt flat)
			: offsets{}, flat(flat)
		{}

		struct equality
		{
			size_t first_false;
			operator bool() { return first_false == Dimensions; }
		};

		[[nodiscard]] constexpr
		equality operator==(const iterator& other) const
		{
			// this is find_if, but pointers make optimizer comfooze
			for(size_t i = 0; i < offsets.dimensions; ++i)
				if(offsets[i] != other.offsets[i])
					return {i};
			if(flat != other.flat)
				return {Dimensions - 1};
			return {Dimensions};
		}

		[[nodiscard]] constexpr
		bool operator!=(const iterator& other) const
		{
			// this is inconsistent with ==, not just because it's an optimization,
			// this is a conjunction, while it should be a disjunction
			// see the corresponding operator in vectirator for explanation
			// now how is this a conjunction?
			// it is the invariant of this iterator that if it reaches another iterator in any dimension it would also reach it in all lower dimensions and x is the lowest, therefore if we're not equals in x we can't be equal in any other dimension. In other words x != is true when all other dimensions != is true as well, making the conjunction true, and if x != is false the conjunction is obviously false, regardless of other dimensions.
			return offsets.x() != other.offsets.x();
		}

		// TODO: relational operators, defining a lexicographic order (can just compare FlatIt), as opposed to partial order of vectors

		[[nodiscard]] constexpr
		decltype(auto) operator*() const { return *flat; }

		constexpr FlatIt operator->() const { return flat; }

		constexpr iterator& operator++()
		{
			++offsets;
			++flat;
			return *this;
		}

		[[nodiscard]] constexpr iterator operator++(int)
		{
			iterator old = *this;
			++(*this);
			return old;
		}

		constexpr iterator& operator--()
		{
			--offsets;
			--flat;
			return *this;
		}

		[[nodiscard]] constexpr iterator operator--(int)
		{
			iterator old = *this;
			--(*this);
			return old;
		}

		[[nodiscard]] constexpr
		difference_type operator-(const iterator& other) const
		{
			return combine_binary_op<scalar_difference_type>(other, std::minus<>{});
		}

		[[nodiscard]] constexpr
		friend iterator operator+(iterator it, const difference_type& diff)
		{
			it.split_binary_op(diff, std::plus<>());
			return it;
		}

		[[nodiscard]] constexpr
		friend iterator operator+(const difference_type& diff, iterator it)
		{ return it + diff; }

		[[nodiscard]] constexpr
		friend iterator operator-(iterator it, const difference_type& diff)
		{
			it.split_binary_op(diff, std::minus<>());
			return it;
		}

		[[nodiscard]] constexpr
		friend iterator operator-(const difference_type& diff, iterator it)
		{ return it - diff; }

		// TODO: in place arithmetic ops

		constexpr void next(const offset_type& step, size_t dimension)
		{
			// dimension == 0 means we are in the beginning or somewhere in the middle, and just need to increment
			assert(dimension > 0);

			// dimension == Dimensions means we are at the end of the range, nowhere to go
			assert(dimension < Dimensions);

			for(size_t i = 0; i < dimension; ++i)
				offsets[i] = 0;

			for(size_t i = dimension; i < offset_type::dimensions; ++i)
				offsets[i] += step[dimension - 1];

			flat += step[dimension - 1];
		}

		constexpr void next(const offset_type& step, equality done)
		{
			next(step, done.first_false);
		}

		constexpr void prev(offset_type step, size_t dimension); // TODO
		constexpr void prev(offset_type step, equality done); // TODO

		constexpr void advance(vector<scalar_difference_type, offset_type::dimensions> step,  size_t dimentins); // TODO
		constexpr void advance(vector<scalar_difference_type, offset_type::dimensions> step, equality done); // TODO

		private:
		offset_type offsets;
		FlatIt flat;

		template <typename Result, typename BinaryOp>
		[[nodiscard]] constexpr
		vector<Result,Dimensions> combine_binary_op(const iterator& other, BinaryOp&& bop) const
		{
			return vector<Result, offset_type::dimensions>(bop(offsets,other.offsets))
				.concat(bop(flat,other.flat));
		}

		template <typename Combined, typename BinaryOp>
		constexpr void split_binary_op(Combined value, BinaryOp&& bop)
		{
			offsets = bop(offsets, value.template first<offset_type::dimensions>());
			flat = bop(flat, get<offset_type::dimensions>(value));
		}

	};

	// a more straight forward version
	// all offset as present, along with the underlying iterator that serves as a reference to the beginning of the range
	// so the logical structure is
	// 2D - vector(row_offset, page_offset), begin_iterator
	// 3D - vector(row_offset, page_offset, book_offset), begin_iterator
	// 4D - vector(row_offset, page_offset, book_offset, shelf_offset), begin_iterator
	// etc
	// advantages:
	// simple implementation (in theory may be easier for compilers to optimize, and any future geom::vector optimizations will adhere)
	// inevitably carries more information, which may be of use in some contexts
	// generic interface
	// disadvantages:
	// slow in practice
	// inevitably carries more information, which makes it more bulky
	// TODO: compatibility with other version geom::iterator
	template <typename BeginIt, size_t Dimensions = 2>
	class vectirator
	{

		public:

		using offset_type = vector<std::uintptr_t, Dimensions>;
		using scalar_difference_type = typename std::iterator_traits<BeginIt>::difference_type;
		using difference_type = vector<scalar_difference_type, Dimensions>;
		using value_type = typename std::iterator_traits<BeginIt>::value_type;
		using reference = typename std::iterator_traits<BeginIt>::reference;
		using pointer = typename std::iterator_traits<BeginIt>::pointer;

		// TODO: constrain BeginIt to not be input or output iterator
		// TODO: use generic iterator function instead of the operators in the implementation
		// NOTE: category for this is always bidirectional or partially random access(see geom::iterator), however if BeginIt is forward iterator dereference will be expensive
		using iterator_category = std::bidirectional_iterator_tag;
		//TODO: iterator_concept

		constexpr vectirator(BeginIt origin, offset_type offsets)
			: offsets(offsets), origin(origin)
		{}
		constexpr vectirator(BeginIt origin)
			: offsets{}, origin(origin)
		{}

		[[nodiscard]] constexpr
		auto operator==(const vectirator& other) const
		{
			return offsets == other.offsets;
		}

		[[nodiscard]] constexpr
		auto operator!=(const vectirator& other) const
		{
			// this is an unfortunate thing with standard loops that they use != and thus in multidimensional case to stop with the shortest sequence we need a conjunction, which is totally counter-intuitive since, when it comes to actually comparing the values, you expect a disjunction
			// as a result we now have this inconsistent with ==, to keep == sane for normal comparison... ideally standard loops conditions would have a customization point other than != that would default to !=, but alas... we are here now, if you want to do logic use !(x==b), x!=b is special customization point reserved for loops... technically i must wreck == too since begin != end is just a common convention not a requirement, but I'll hold on to it for as long as I can
			// why do we want to support standard loops anyway? because in our case they correspond to a loop over a row, which is a very common thing, and with large arrays is the part that actually matters in terms of performance, so that's where all the fancy algorithms will come in while the outer higher dimensional loops can be their own thing using next(step, mask) et al
			// TODO: on the other hand since there is always going to be a fancy outer loop, we can switch types for the inner loop and keep == and != consistent, sane for the outer loop  and outside world in general and insane only for the short lived inner loop
			return to_conjunction(offsets != other.offsets);
		}

		// TODO: relational operators (as geom::iterator)

		constexpr BeginIt operator->() const { return origin + *(offsets.end()-1); }

		[[nodiscard]] constexpr
		decltype(auto) operator*() const { return *(operator->()); }

		constexpr vectirator& operator++()
		{
			++offsets;
			return *this;
		}

		[[nodiscard]] constexpr vectirator operator++(int)
		{
			vectirator old = *this;
			++(*this);
			return old;
		}

		constexpr vectirator& operator--()
		{
			--offsets;
			return *this;
		}

		[[nodiscard]] constexpr vectirator operator--(int)
		{
			vectirator old = *this;
			--(*this);
			return old;
		}

		[[nodiscard]] constexpr
		difference_type operator-(const vectirator& other) const
		{
			return offsets - other.offsets;
		}

		[[nodiscard]] constexpr friend
		vectirator operator+(vectirator it, const difference_type& diff)
		{
			it.offsets += diff;
			return it;
		}

		[[nodiscard]] constexpr friend
		vectirator operator+(const difference_type& diff, vectirator it)
		{ return it + diff; }

		[[nodiscard]] constexpr friend
		vectirator operator-(vectirator it, const difference_type& diff)
		{
			it.offsets -= diff;
			return it;
		}

		[[nodiscard]] constexpr friend
		vectirator operator-(const difference_type& diff, vectirator it)
		{ return it - diff; }

		constexpr vectirator& operator+=(const difference_type& diff)
		{ offsets += diff; return *this; }
		constexpr vectirator& operator-=(const difference_type& diff)
		{ offsets -= diff; return *this; }

		constexpr void next(const offset_type& step, vector<bool, Dimensions> done_mask)
		{
			auto begin = std::begin(done_mask);
			auto end = std::end(done_mask);
			auto first_false = std::find(begin, end, false) - begin;
			assert(first_false < step.dimensions);

			offsets += step[first_false];
			offsets *= ~done_mask;
		}

		constexpr void prev(offset_type step, vector<bool, Dimensions> done_mask); // TODO

		constexpr void advance(const offset_type& step, vector<bool, Dimensions> done_mask); // TODO

		private:
		offset_type offsets;
		BeginIt origin;
	};

} // namespace simple::geom

#endif /* end of include guard */
